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VO ! Plateaux in the magnetization curves of the square, triangular and hexagonal lattice spin- 

1/2 XXZ antiferromagnet are investigated. One finds a zero magnetization plateau (cor- 
^ ■ responding to a spin-gap) on the square and hexagonal lattice with Ising-like anisotropies, 

0^ , and a plateau with one third of the saturation magnetization on the triangular lattice 

which survives a small amount of easy-plane anisotropy. Here we start with transfer ma- 
c^ I trix computations for the Ising limit and continue with series in the XXZ-anisotropy for 

C ' plateau-boundaries using the groundstates of the Ising limit. The main focus is then a 

i-^ ■ numerical computation of the magnetization curves with anisotropies in the vicinity of the 

C , isotropic situation. Finally, we discuss the universality class associated to the asymptotic 

behaviour of the magnetization curve close to saturation, as observed numerically in two 
and higher dimensions. 



A Feodor-Lynen fellow of the Alexander von Humboldt-foundation. 



1. Introduction 



The discovery of high-Tc superconductivity has revived interest in two-dimensional Heisen- 
berg antiferromagnets (for a review see [1]), since the Cu02 planes give rise to a good 
realization of the S = 1/2 square lattice antiferromagnet. Due to the large coupling con- 
stants of the high-Tc materials, the main focus is on properties in zero or small external 
magnetic fields. Nevertheless, it has also been progressively realized that antiferromagnets 
exhibit interesting phenomena in strong external magnetic fields, namely plateaux in their 
magnetization curves at certain fractions of the saturation value of the magnetization. 

In one dimension, the appearance of plateaux in magnetization curves is by now rather 
well understood in terms of a quantization condition on the magnetization that involves 
the volume of a translationally invariant unit cell [2 — 8] . Actually, if the interaction inside 
finite clusters of spins is large with respect to the other interactions, the appearance of 
plateaux is governed by the volume of such a cluster irrespective of the dimension [8]. 
However, the simplest two-dimensional systems have equal coupling constants, and then 
it is less clear what determines the appearance of plateaux in magnetization curves. One 
well-known example of such plateaux in two dimensions is a plateau at one third of the 
saturation magnetization in the triangular lattice antiferromagnet. There is in fact a 
number of theoretical studies of the magnetization process of two-dimensional triangular 
antiferromagnets (see [9 — 15] for a selection) which are at least to some extent motivated 
by the presence of this plateau or the more general feature of frustration. The number 
of recent investigations of the magnetization process of the square lattice antiferromagnet 
[16 — 20] still seems to be smaller, and for the hexagonal lattice we are aware of just a 
single study of the Ising antiferromagnet in a magnetic field [21]. 

The magnetization plateau of the triangular lattice antiferromagnet can be observed ex- 
perimentally in a variety of materials (see e.g. [22] for recent rather clear examples and 
[23] for a review of experimental facts about triangular lattice antiferromagnets). The 
present investigation was in fact originally motivated by a magnetization experiment on 
the stacked triangular lattice antiferromagnet CsCuCls [24] which shows a plateau in the 
magnetization curve at one third of the saturation value if the field is applied perpendic- 
ular to the c-axis. In many cases and in particular in CsCuCls, the spin is carried by a 
Cu^+ ion, giving rise to a spin S = 1/2. Furthermore, at least in CsCuCls, anisotropy of 
the interaction is important. This lead us to considering a spin-1/2 XXZ model and to 
investigate the effect of the XXZ-anisotropy A. The relevance of the results of the present 
paper for CsCuCla will be discussed elsewhere [25]. 

The focus of the present work are general features, namely in which situations magne- 
tization plateaux arise in two dimensions and the universality class of the transition to 
saturation. These questions will be investigated by computing the zero-temperature mag- 
netization process of the 5" — 1/2 XXZ model on the aforementioned three lattice types, i.e. 
on a square, triangular and hexagonal lattice. This problem is described by the following 
Hamiltonian: 

H = JJ2 { A^l^l + h {siS^ + SZS+) }-hJ2S^,, (1.1) 

{x,y) X 



where the 5"^ are spin-1/2 operators acting at place x and /i is a dimensionless magnetic 
field. The notation (x, y) denotes neighbouring pairs on a lattice whose total number of 
sites we denote by V. The magnetization (M) is given by the expectation value of the 
operator M = y ^^ S^ where the prefactor is chosen in order to normalize the saturation 
value to (M) = ±1. The magnetization operator M commutes with the Hamiltonian (1.1). 
This leads to a technically useful simplification since it allows one to relate all properties 
in a magnetic field h to those at h = with a suitably fixed magnetization (M). 

The plan of this paper is as follows: First, we compute magnetization curves for the Ising 
model which is obtained from (1.1) by dropping the S^S~ hopping matrix elements. Exact 
zero-temperature groundstates are readily written down for all plateaux observed in this 
limit. We then use this as an input to compute perturbation series in A~^ for the gap of 
single-spin excitations above these groundstates. Even though the boundaries of plateaux 
are in general determined by multi-spin excitations, the series for the single-spin excitations 
yield a guide in which region the plateaux persist for A < oo. For general XXZ-anisotropy 
A, plateaux existing in the Ising limit A = oo may not just disappear, but further plateaux 
could arise. This and the intrinsic limitations of a perturbative approach necessitates a 
direct computation of the magnetization process of the full quantum Hamiltonian (1.1), 
in particular in the region where A is of order unity. The bulk of the paper is therefore 
devoted to a numerical investigation of (1.1) on all three lattice types which is much in 
the spirit of classical work on single Heisenberg chains [26,27]. Finally we discuss the 
universality class associated to the asymptotic behaviour of the magnetization curve close 
to saturation, i.e. for (M) -^ 1. 

2. The Ising antiferromagnet 

The simplest case where one can observe magnetization plateaux in two dimensions is the 
Ising antiferromagnet. The latter can be obtained from (1.1) by taking the limit A ^ oo 
and rescaling J -^ J/ A. This yields the energy function 

{x,y) X 

where the sg are now Ising variables with s^ = ±1. Note that here we use unusual con- 
ventions for J and h in order to later simplify making contact with the XXZ Hamiltonian. 

While in the bulk of the paper we restrict to zero temperature, here it is actually useful to 
work at finite inverse temperature P in order to permit application of the transfer-matrix 
method (see e.g. [28]). The magnetization is then given by 

iM) = ^' :-,.„..„ ■ (2-2) 

Note that here M is no longer a conserved quantity and one can therefore not drop the 
expectation values. Working with expectation values would also be mandatory for an 
anisotropy axis which does not coincide with the field direction, but this is not considered 
in the present paper. 



In principle, one could directly write down zero-temperature groundstates of (2.1) and 
compute the values of h where one of them becomes preferred to others. Actually, a 
complete set of groundstates is known since the two-dimensional Ising model in the presence 
of an external magnetic field has already been studied some time ago on the square lattice 
(see e.g. [29]), the triangular lattice (see e.g. [30]) and somewhat more recently on the 
hexagonal lattice [21]. 

However, since we later wish to build on the results for this limit, we believe that it is still 
useful to summarize the results in a uniform framework. We therefore evaluate (2.2) on a 
strip using the transfer-matrix method [28], where we employ periodic boundary conditions 
along the the short direction of the strip and open ones in the long one. The site y in (2.2) 
is put at the center of the strip. 
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Fig. 1: Magnetization curves of an Ising antiferromagnet on a 18 x 1000 
a) square, b) triangular and c) hexagonal lattice. The lines are for 

J/3 = oo (full), J/3 = 40 (long dashes) and J {3 = 8 (short dashes). 



Fig. 1 shows magnetization curves obtained in this manner on a square, triangular and 
hexagonal lattice. The strip for Fig. 1 was always chosen to be 1000 sites long and 18 sites 
wide. A computation with the same length but half the width (i.e. 10 sites for the square 
and hexagonal lattice, but 9 for the triangular one) leads to values for (M) which differ 
at most by 10~^ from the ones shown. We may therefore expect that the magnetization 
curves in Fig. 1 are basically indistinguishable from those of the thermodynamic limit. 
The smallness of finite-size effects will be a justification for using small system sizes later 
in the analysis of the XXZ antiferromagnet, although then finite-size corrections should be 
considered. 

The interpretation of Fig. 1 is as follows. Both the square and hexagonal lattices are 
bipartite lattices. Thus, their zero-temperature groundstate at sufficiently small fields is 
given by all spins pointing up on one sublattice and pointing down on the other. If the 
magnetic field becomes large enough to make it favourable to align one further spin along 
the field, so it becomes for all and one has a sharp transition from the unmagnetized 
antiferromagnetic groundstate to a fully magnetized ferromagnetic one. The transition 
field is readily computed to equal /i = 2J or /i = | J for the square and hexagonal lattice. 



respectively (actually, it is easy to see that the transition to a fully magnetized state takes 
place at /i = I J, where z is the coordination number of the lattice). 

The situation is slightly different for the triangular lattice where each plaquette is frus- 
trated. This lattice has three sublattices and energy is minimized at small fields by aligning 
all spins in two of them along the field and the ones in the third in the opposite direction. 
Increasing the field it becomes again favourable for all the latter spins simultaneously to 
align also along the field. Thus one obtains a sharp first-order transition from a state with 
(M) = 1/3 to the fully magnetized one at /i = 3J. 

These transitions are sharp only at zero temperature (/? = oo) and smoothed out by thermal 
fiuctuations otherwise. Since such a smoothing effect of finite temperature is generic, we 
will restrict from now on to zero temperature after having illustrated this effect in the case 
of the Ising model. 

3. Expansions around the Ising limit 

After this illustrative study of the Ising antiferromagent we now turn to the XXZ model 
(1.1). We use the groundstates described in the previous section to expand the gap of 
single-spin excitations in powers of A~^. 

Since we wish to cover a variety of cases, it is convenient to use a simple but general method 
for higher order series expansions of a quantum mechanical system which is summarized 
e.g. in Section 3 of [31]. This should be sufficient to obtain an overview, but could certainly 
be extended to higher orders using more sophisticated cluster expansions if this should turn 
out to be desireable for concrete applications. 

For the square and hexagonal lattice there is just a plateau at (M) = 0. The lowest 
excitations are those where a single spin is fiipped with respect to the antiferromagnetic 
groundstate. Due to the antiferromagnetic nature of the groundstate, the first-order cor- 
rections to the energies vanish and one finds a non-trivial dispersion only in second order 
in A~^. Since both cases have been studied in detail in [32] (see also references therein) 
and [33] , respectively, we skip the details of the computation. 

Analytical fourth-order expressions for this 5"^ = 1 gap are given by 

^ = -^ f 2A - A + J|L^ + o(A-4) (3.1) 



3A 432A3 



for the square lattice, and 



/3, 15 295 \ ,^,, 4, , , 

for the hexagonal lattice. Numerical versions of the coefficients up to tenth order can be 
found in Table I of [32] and [33], respectively. Applying a naive ratio test to the higher 
orders given in [32], one concludes that the series (3.1) can be expected to converge for 
A-2<1. 



Now we proceed with the triangular lattice. RecaU that this has three sublattices and the 
lowest-energy state with (M) = 1/3 is obtained by ahgning aU spins on two sublattices up 
and the ones on the third one down. The energy of this state at /i = is found to be 

One possible excitation is obtained by flipping one further spin up. Using a Fourier trans- 
formation to lift the translational degeneracy, one finds the dispersion relation 

r cos kj. + cos ky — ^ — 2 (cos k^ + cos ky) cos(kjr + ky) — cos(kjr — ky) 
S+{k^,ky)=J 3A+ ^ ^ ^ — ^- — ^ — ^ 

+ C(A-2). 

(3.4) 
From its minimum one obtains the gap for dS^ = 1 excitations above the (M) = 1/3 
plateau 

+ ^ ' ^ + V 3 3 y 4A 16A2 320A3 ^ ' ^ ' 

Now consider excitations with S^ by 1 smaller than on the (M) = 1/3 plateau state. There 
are two possibilities to fiip a spin down with respect to the plateau state which are not 
related by translational symmetry. In order to minimize the energy one has to take the 
difference as a linear combination of these two possibilities. Then one finds in a way similar 

to (3.4) 



^— (fex7 ky) J 



i(cosfc. + cosfc, + cos(fc. + M) 

1 /7 

-|- — I — h cos(A;a; + A;y)(cos k^ + cos ky — cos k^ cos ky) 
A y 8 

3 

— -{cos{kx + ky) + cos kx + cos ky) + cos k^ cos ky 



+ C(A-2). 



(3.6) 
Its minimum determines the gap for 5S^ = — 1 excitations above the (M) = 1/3 plateau 
and is given by 

^ . X 3J 5J 73J 42787J ,^, , ., , „, 

These series will be compared to results of a numerical diagonalization in the following 
sections. 

We conclude this section by mentioning that the upper critical field hyc at which the 
transition to a fully magnetized state takes place is straightforwardly computed exactly if 
it is determined by a single spin-fiip. One finds 

{dJ{A + 1) ((i-dimensional hypercubic lattice), 
3 J(A + |) (triangular lattice), (3-8) 

|j(A + 1) (hexagonal lattice). 



For later use we have given here actuaUy the value of hue for the (i-dimensional hypercubic 
lattice - the result for the square lattice is given by its (i = 2 special case. Note that 
the overall numerical factor is proportional to the coordination number z of the lattice. 
Furthermore, for the bipartite lattices (hypercubic and hexagonal) one has hue = ^J{^ + 

4. Numerical diagonalization for the square lattice 

In this and the following sections we report results of a numerical study of the magneti- 
zation process for XXZ anisotropies around A = 1. We have numerically calculated the 
lowest eigenvalues as a function of the magnetization, wave vectors and A on finite sys- 
tems. Below we just present the consequences for the magnetization curves and skip the 
details of these standard (but still CPU time intensive) computations. 

We choose to present our results in terms of 'magnetic phase diagrams'. They show the 
projection of the more conventional magnetization curves onto the axis of the magnetic 
field. The values of (M) are thus assigned to different regions of the plot. Since in this 
compact representation we save one axis, we can display the variation of the magnetization 
curves as a function of A in a single figure. 

First we discuss the square lattice. This case has also been studied with finite system 
diagonalizations atA = lona4x4 lattice [16] as well as more recently on larger lattices 
[17]. At A = 1 also a second-order spin- wave investigation has been performed [19]. For 
A > 1 finite system diagonalizations and quantum Monte Carlo simulations have been 
carried out in [18]. The system is therefore well understood and provides a good check 
of our method. Before presenting our results we recall from [18] that one finds a plateau 
with zero magnetization for A > 1 whose boundary corresponds to a first-order phase 
transition, i.e. in the thermodynamic limit the magnetization jumps by a finite amount. 
This plateau disappears (its width tends to zero) as A ^ 1. 

Fig. 2 shows the magnetic phase diagram in the region with A close to one on a square 
lattice of size 4x6, i.e. with a volume of 24 spins. The thin full lines denote boundaries 
of the magnetization plateaux (M) = to/12 (to = 0, . . . , 12) which have to occur for this 
system size. For other system sizes also other values for (M) will be possible. Therefore, 
regions in Fig. 2 where these lines are regularly spaced can be expected to correspond to 
smooth transitions in the thermodynamic limit. Bearing this in mind, one can clearly see a 
plateau with magnetization (M) = for A > 1.05. On the other hand, this plateau seems 

to be absent for A < 0.95. It should be noted that for A > 1.325 the neighbour of (M) = 

has magnetization 1/6, i.e. here it is favourable to fiip two spins in the direction of the 
magnetic field rather than one. This fact refiects the first-order nature of the transition. 
In this region A > 1.325, the energy required for an excitation corresponding to a single 
fiipped spin is shown by a dashed line (otherwise the energy for such an excitation is equal 
to the boundary of the (M) = plateau: hie = "^)- 

The extension to tenth order [32] of the series (3.1) is shown by the the bold full line in 
Fig. 2. It should be compared to the energy-gap for a single fiipped spin which is not 
always equal to the boundary of the (M) = plateau (in this case shown by the dashed 

6 



line). In the interval 1.1 < A < 1.5, this series slightly overshoots the finite-size data. In 
fact, the same can be observed in earlier presentations [32, 18] and therefore probably is 
not due to a finite-size effect but missing higher-order terms in the series which could be 
quite important since we are looking at a region close to the limits of validity of this series 
( A > 1 should be its region of convergence) . 
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Fig. 2: Magnetic phase diagram of the square lattice XXZ antiferro- 
magnet on a 4 x 6 lattice. The thin full lines denote boundaries of areas 
with the values of the magnetization indicated in the figure. The dashed 
line is the spin-gap (see text). The bold full line shows the extension to 
tenth order [32] of the series (3.1) for the gap. 



The fully magnetized state with (M) = 1 gives rise to a further trivial plateau in Fig. 2. 
The finite-size data for its boundary agrees with the analytical result (3.8) as it should. 

Let us summarize our results for the square lattice: Using a system of size 4 x 6 we are able 
to locate the end of the zero magnetization plateau in the interval 0.95 ^ A < 1.05 which 

is a reasonably good approximation to the presumably exact value A = 1. Furthermore, 
the first-order nature of the transition refiects in the fact that with this lattice size it is 
favourable to fiip two spins rather than one for A > 1.325. 



5. The triangular lattice numerically 

The triangular lattice was already studied some time ago by finite system diagonalizations 
[12]. Now we can access larger system sizes and look more carefully at the dependence on 
A ^). To obtain an overview, we first look at the 3x6 lattice, even though precisely this 
lattice was already studied in [12] for A G {0.8, 1, 1.2, 2.5, 5}. Fig. 3 shows our magnetic 
phase diagram for this case. The magnetization curves at A = 1, 1.2 and 2.5 of [12] 
correspond to sections through Fig. 3 and are in all cases in excellent agreement with our 
results. 
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Fig. 3: Magnetic phase diagram of the triangular lattice XXZ antifer- 
romagnet on a 3 x 6 lattice. The thin full lines denote boundaries of 
areas with the values of the magnetization indicated in the figure. The 
dashed line is a single-spin excitation above the (M) = 1/3 groundstate 
(see text). The bold full lines show the series (3.5) and (3.7). 



Here we observe precisely one non-trivial plateau with (M) = 1/3. Note that for A > 2.4 it 
becomes favourable to fiip two spins rather than one at the upper boundary of this plateau. 



^) The current record for finite system diagonalizations on this type of lattice seems to be 
held with a volume of 36 spins (see e.g. [34, 14]). Since we want to vary the magnetization, 
wave vectors, A and still limit the computational effort, we content ourselves with smaller 
system sizes. 



In this case, the energy of an excitation corresponding to a single spin flipped is shown 
by a dashed line while the finite-size data for transitions between different groundstates is 
shown by thin full lines. Based on our experience with the square lattice we take this as 
an indication that the transition at the upper boundary of this plateau becomes first order 
for A > 2 while all other transitions are second order. 

The series (3.5) and (3.7) for the upper hc^ = ^+(0,0) and lower boundaries hc^ = 
— £^_(0, 0) of the (M) = 1/3 plateau are shown by bold full lines. The former should 
be compared to the dashed line, the latter to the appropriate thin full line. The agreement 
is good for the right half of Fig. 3. Since here we have less orders for the series than for the 
square lattice (Fig. 2), it is not suprising that the agreement in the region A close to one 
is less good. Again, the finite-size data for the location of the transition (M) — ^ 1 agrees 
exactly with (3.8). 
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Fig. 4: Magnetization curve of the spin-1/2 triangular lattice Heisen- 
berg antiferromagnet with A = 1 on 3 x 6 (short dashes), y = 27 (long 
dashes) and V = 36 (thin full line) lattices. The curves with volume 
y = 27 and y = 36 for (M) < 1/3 are based on data of [14]. The bold 
full line is an extrapolation to the thermodynamic limit (see text). 

Now we examine the region A close to one in more detail. First, we present the magnetiza- 
tion curve at A = 1 in Fig. 4. The thin lines denote curves at three different system sizes. 
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Here, we have used results of [14] to obtain the parts with (M) < 1/3 of the magnetization 
curves with 1/ = 27 and V = 36. Our results overlap with those of [14] just at the lower 
boundary of the (M) = 1/3 plateau at a volume V = 27. However, our geometry is differ- 
ent from that of [14]: In our computation the V^ = 27 lattice has a smallest spatial extent of 
just 3 sites while the configuration of [14] was designed to maximize the distance between 
boundaries. This leads to a difference of 2% between the results, which is reasonably small 
in view of the small linear size. 

In Fig. 4 one can clearly see a plateau with (M) = 1/3. The finite-size effects for its 
boundaries are small ^). One observes that finite-size effects are also small for the midpoints 
of the steps in these curves, as is well-known from one dimension [26,27]. Drawing a 
curve through these points for the largest available system size (here we have used up to 
V = 225 = 15 X 15 in the vicinity of the upper critical field hue), we therefore obtain the 
approximation to the thermodynamic limit shown by the bold full line in Fig. 4. 
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Fig. 5: Magnetic phase diagram of the triangular lattice XXZ antifer- 
romagnet on a 3 x 6 (full line), 3x9 (long dashed line) and 6x6 (short 
dashed line) lattice. 



^) The first-order spin- wave results for the boundaries of the (M) = 1/3 plateau, /ici ~ 
1.248J and hc^ ~ 2.145J [13] are about 0.13J smaller than the finite-size diagonalization 
results. 
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Finally, in Fig. 5 we zoom in to the region A around one of Fig. 3, using also bigger system 
sizes "^) (Fig. 4 is a section at A = 1 through Fig. 5). The main motivation for taking a 
closer look at this region comes from the observation in [12] that the (M) = 1/3 plateau 
is present at A = 1, but does not seem to exist at A = 0.8. Inspecting Fig. 5 and paying 
attention to the finite-size effects (in particular the difference between the 3x6 and the 
3x9 lattice), one concludes that this plateau presumably disappears somewhere in the 
region A ^ 0.85. 

We have also indicated the location of the (M) = 2/3 plateau which is only realized for 
the 3x6 and the 6x6 lattice. However, the finite-size data provides no indication for it 
to survive in the thermodynamic limit. So, it (and possible other plateaux) are likely to 
be absent in the thermodynamic limit, as is implied by the bold full line in Fig. 4. 

In Fig. 5 we have omitted the series (3.5) and (3.7) since already in Fig. 3 one can observe 
notable deviations in this region of A. Of course, the result (3.8) is still valid. 

6. Numerical diagonalization for the hexagonal lattice 

While the square lattice is self-dual, the dual of the triangular lattice is the hexagonal 
lattice. It is therefore interesting to also investigate the magnetization process on the 
hexagonal lattice. Fig. 6 shows the results of diagonalizations on a 4 x 6 lattice. The 
region of small magnetization looks qualitatively very similar to that of the square lattice 
in Fig. 2. In particular, one can see an (M) = plateau (corresponding to the gap) for 
A>1. Inferring from Fig. 6 that there is a first-order phase transition at the boundary 

of this plateau may be speculative, but this should also be expected on the grounds of 
universality, i.e. in the thermodynamic limit the transition at the boundary of the (M) = 
plateau for the hexagonal lattice should be in the same universality class as that of the 
square lattice. The location of the ending-point of the (M) = plateau is compatible with 
A = 1 (see also [33]) which would be the same as that for the square lattice. While in 
general one cannot use arguments based on universality to locate a critical point, here the 
point A = 1 is distinguished by enhanced su{2) symmetry and one may therefore expect 
the closing of the gap exactly at A = 1. 

The tenth order version of [33] of the series (3.2) is shown by the bold full line in Fig. 6 
and should be compared with the dashed line. The agreement is quite good close to the 
right boundary of the figure and less good for smaller values of A. This is not surprising 
since the series for the hexagonal lattice [33] clearly converges less well than that for the 
square lattice [32]. 

As far as non-zero magnetizations are concerned, observe first that the transition to the 
fully magnetized state does indeed take place at the value of hue given by (3.8), thus giving 
a crosscheck on our computations. More important are the candidates for plateaux. In 
Fig. 6, the only plausible value for the appearance of a non-trivial plateau is at (M) = 1/2. 



"^) We have not computed all transition lines for the larger lattices sizes. We have 
omitted the ones for (M) < 1/3 on the 3x9 lattice and those with (M) < 2/3 on the 6x6 
lattice. 
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Fig. 6: Same as Fig. 2, but for ttie tiexagonal lattice. The bold full line 
shows the extension to tenth order [33] of the series (3.2). 



It should be noted that the geometry corresponding to Fig. 6 can be interpreted as a 
variant of the A^ = 4-leg spin-ladder, and here plateaux are expected at (M) = and 
(M) = 1/2 [4,5]. So, in order to see whether the wide step in Fig, 6 at (M) = 1/2 is 
an intrinsic feature of the hexagonal lattice, one should look at other system sizes, and 
in particular larger strip widths. Since the large computational effort makes a systematic 
investigation of the dependence on A unfeasible, we have looked at the 6x6 lattice for 
a few selected values of A. The finite-size magnetization curves for A = 1 in Fig. 7 are 
representative of the general case. In addition to the data for the 4x6 lattice, we here 
also show magnetization curves for the 4x4 lattice and 6x6 lattice (close to saturation, 
some larger lattice sizes have also been considered, but are not explicitly shown). 

With the data for the 6x6 lattice taken into account, there is no indication anymore of 
an (M) = 1/2 plateau. The bold full line therefore shows an extrapolated magnetization 
curve without plateaux, which as in Fig. 4 was obtained by connecting the midpoints of 
the steps for the largest available system sizes. This extrapolated curve shows that the 
4x6 lattice is still subject to substantial finite-size effects which at (M) = 1/2 conspire 
to suggest a plateau. 

Also at A = 2 there is no evidence either for a plateau with (M) = 1/2 on a 6 x 6 lattice 
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(note that here the corresponding candidate in Fig. 6 on a 4 x 6 lattice is most pronounced). 
Thus, the final picture for the hexagonal lattice is the same as for the square lattice, i.e. 
the only plateau occurs for A > 1 at (M) = 0. The main qualitative difference between 
these two lattice types is in the finite-size corrections which are much more important for 
the hexagonal lattice. 
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Fig. 7: Magnetization curve of the hexagonal antiferromagnet with 
A = 1. The finite-size data is for 6x6 (thin full line), 4x6 (long 
dashes) and 4x4 (short dashes) lattices, respectively. As in Fig. 4, the 
bold full line is an extrapolation to infinite volume. 



7. The transition to saturation 



Finally we look at the asymptotic behaviour the magnetization as a function of the mag- 
netic field close to saturation. One possibility is that it reaches the upper critical field hue 
with a power law. Then one would introduce a critical exponent n via 



1-(M) 



h. 



h 



J 



A different possibility is a linear behaviour with a logarithmic correction 
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(7.2) 
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The latter has been argued in [19] to apply to the square lattice antiferromagnet (see also 

[35,36]). 

To determine the functional form and estimate n we use the midpoint of the last step in 
the magnetization curve, as was done for the one-dimensional case e.g. in [37]. Since the 
corresponding magnetization satisfies 1 — (M) = 2/V, determination of the asymptotic 
behaviour of the magnetization curve is equivalent to determining the asymptotic finite- 
size behaviour of this last step. The particular choice of the midpoint of the step is not 
relevant for the determination of the exponent (or more generally, the functional form). 
Any other choice such as the last corner would simply yield a different prefactor. Values of 
such midpoints are given in Table 1 for all three lattice types at the isotropic point A = 1. 

Fitting this data for the square lattice to the form (7.1), we find a value of n in the region 
H ~ 0.83. However, the precise value increases if we use larger system sizes for the fit, as 
is expected to be the case in the presence of a logarithmic correction. Indeed, one obtains 
a better fit if one uses (7.2) instead. Using all data points for the square lattice in Table 
1 we find a = 0.2505 ± 0.0015 and b = 2.148 ± 0.075 where the errors indicate the la 
confidence interval of the fit. So, the numerical data can be regarded as a confirmation of 
the functional form (7.2) predicted by first-order spin- wave theory [19], though it is not 
surprising that the values for the constants differ from the first-order spin-wave predictions 
[19] which specialized to 5" = 1/2 are a = l/(27r), b = tt^. In no case do we reproduce 
the simple linear behaviour reported in [17]. The crucial difference is probably not that 
we employ exclusively system sizes which are larger than those used in [17], but that the 
analysis in this reference not only assumed the form (7.1), but also that 1/// is an integer. 



square lattice 


triangular lattice 
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0.028767 
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0.057168 
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0.022087 
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0.0089914 
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0.010450 


12 X 12 


0.0094961 


12 X 12 


0.0057555 


14 X 14 


0.0071704 


15 X 15 


0.0054802 


14 X 14 


0.0039638 


16 X 16 


0.0051911 


18 X 18 


0.0035220 


16 X 16 


0.0028778 


18 X 18 


0.0039129 


21 x21 


0.0024338 


18 X 18 


0.0021741 


20 X 20 


0.0030438 


24x24 


0.0017719 


20 X 20 


0.0016943 


22 X 22 


0.0024282 


27 X 27 


0.0013419 


22 X 22 


0.0013540 


24 X 24 


0.0019776 


30x30 


0.0010480 


24 X 24 


0.0011045 


26 X 26 


0.0016386 


33x33 


0.00083883 


26 X 26 


0.00091571 



Table 1: Midpoint of the last step before the upper critical field hue 
for two-dimensional antiferromagnets at A = 1. 

In the case of the triangular lattice, one should first discard the data for he volumes 
y = 6 X 6 and V = 9 x 9 in order to obtain a smooth line. Then one again obtains the 
same behaviour as for the square lattice: The exponent /i obtained by the fit (7.1) has a 
very similar value and moves in the same direction with increasing V as for the square 
lattice. Furthermore, one obtains a better fit with (7.2) than with (7.1). Thus we conclude 
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that also the transition to saturation of the triangular lattice antiferromagnet obeys (7.2), 
where the constants are now given by a = 0.1495 ± 0.0006, b = 2.498 ± 0.061. 

For the hexagonal lattice we can access the same system sizes as for the square lattice. 
Thus, we can directly compare the values for the midpoints of the steps. One finds that 
those for the hexagonal lattice differ from the ones for the square lattice by a factor which 
approaches rapidly ~ 0.56 with increasing system size. This means that the transitions 
to saturation of the square and hexagonal lattice should belong to the same universality 
class. Indeed, if we perform the same analysis as for the square lattice, we find the same 
behaviour for the hexagonal lattice. Just the constants for the fit (7.2) are different: Now 
we have a = 0.4338 ± 0.0024, b = 1.527 ± 0.052. 

In summary, we find support for the asymptotic behaviour (7.2) in all cases we have 
considered. This universal behaviour can be understood in a way very similar to the 
DN-PT universal square root in one dimension [38, 39]. One starts from the single-particle 
excitations. At the transition field one starts to fill the lowest band of magnetic excitations. 
Generically, the dispersion around such a minimum is quadratic, e.g. E ~ |/cp. Then one 
needs to know how many states are available below this value of |A;|. In one dimension, 
where the detailed nature of the excitations does not matter, this number is proportional 
to k as long as one has an exclusion principle. Equivalently, in d = 1 the number of states 
available below a given energy E is proportional to y/E. Since the magnetic field h acts as 
a chemical potential and the number of particles corresponds to the deviation of (M) from 
its critical value, this argument leads to the DN-PT universality class in one dimension. 

For hard-core bosons in two dimensions (which is the situation which we consider here), 
interactions lead to a logarithmic correction to the naive dimensional analysis [35,36,40] 
(see also Chapter 6 of [41]) and thus to (7.2). Note however that this argument is indepen- 
dent of the details of the model under consideration. The crucial ingredients are just that 
the fundamental excitations are bosons with a quadratic dispersion around the minimum 
subject to a repulsive interaction. 

In three (and higher) dimensions we have Bose condensation. Therefore, the magnetization 
curve of a hypercubic antiferromagnet in d > 3 should have a simple linear approach to 
saturation [42]. Since this is the behaviour found for classical spins, the mean-field result 
for |U is in some sense exact in d > 3. The effect of Bose condenstation which lowers ^ 
substantially below the value d/2 can be observed clearly in numerical diagonalization of 
the isotropic (hyper)cubic antiferromagnet (compare Table 2). In fact, if we fit the data 
in Table 2 to the form (7.1), we find ^ ~ 0.95 in d = 3 and // ~ 0.98 in d = 4. Given 
the small linear extent of the systems considered, this is in reasonable agreement with the 
predicted value n — 1 without logarithmic corrections. 

Different values for /i can be obtained if the form of the dispersion close to the minimum is 
not quadratic. For example, for special values of parameters one could have i? ~ |/c|^ which 
in d = 1 leads to /x = 1/4. The value fx = 1/4 has been observed in the two- leg zig-zag 
ladder [43] and the biquadratic spin-1 chain [44] with each a special value of the coupling 
constant. This is indeed explained by the aforementioned change in the single-particle 
dispersion (see [45] and [46], respectively). 
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0.00078466 






183 


0.00054581 







Table 2: Midpoint of the last step before the transition to saturation 
for the cubic and d = 4 hypercubic antiferromagnet at A = 1. 

8. Discussion and conclusions 

In the present paper we have observed the foUowing plateaux in magnetization curves: For 
the square and the hexagonal lattice, i.e. the two bipartite lattices, we find an (M) = 
plateau for A > 1. The transition at the boundary of this plateau is likely to be always 
of first order. On the triangular lattice one finds a plateau with (M) = 1/3 for A > 0.85. 

The transitions at its boundary appear to be second order for small enough anisotropies 
(at least for A<2). 

Recently plateaux with (M) = 1/2 (and also (M) = 0) have been observed on the trian- 
gular lattice with multi-spin interaction, in particular a four-spin interaction in addition 
to the two-spin interaction discussed in the present paper [47,48]. At least the (M) = 1/2 
plateau even survives the classical limit [49] 4). This is to be contrasted with the (M) = 1/3 
plateau in the triangular lattice antiferromagnet which is absent in the classical limit and 
arises only in first-order spin- wave theory [13]. 

While all the aforementioned plateaux with (M) 7^ occur in frustrated systems, frus- 
tration is certainly not a necessary ingredient for the appearance of non-trivial plateaux. 
Consider for example a spin-1/2 A^-layer square lattice Heisenberg antiferromagnet (see 
[20] for a detailed discussion of the bilayer system in a magnetic field) . In the limit where 
the in-plane coupling tends to zero, this system decouples into clusters of A^ spins. Simply 
by counting the possible states of these clusters, we conclude that plateaux exist with 
(M) = —1, —1 + 2/A^, . . ., 1 — 2/A^, 1 at least if the inter-plane coupling is much larger 
than the in-plane coupling [4, 5, 8]. However, frustration appears to favour the appearance 
of further plateaux. Therefore, a bilayer triangular lattice may be an interesting candidate 
for further study. 

The triangular and hexagonal lattice are dual to each other and thus share the same 
(local) point-symmetry. However, we have seen that an antiferromagnet in a magnetic 



4) In addition, also a plateau with (M) 
a suitable choice of parameters [49] . 



1/3 can be observed in the classical model for 
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field behaves very difi^erently on them. The point-symmetry group therefore does not 
appear to be of any relevance to the magnetization process. On the other hand, the square 
and hexagonal lattice share the property of being bipartite lattices and in fact gives rise 
to very similar behaviour in the presence of a magnetic field. 

In the classical or Ising limit it is clear that the appearance of plateaux is related to the 
number of sublattices needed to describe the full magnetization process. It is therefore 
likely that some kind of translationally invariant unit cell will control the appearance of 
plateaux also in higher dimensions. The definition of such a unit cell in two dimensions is 
nevertheless far less obvious than in one dimension, and some ambiguity is also possible. 
To elucidate the situation further, it would be useful to study a wider range of models and 
lattice types. An important class of lattice types would be given by the eleven Archimedean 
tilings (see e.g. Chapter 2 of [50]) among which we have considered the three monohedral 
ones. In particular, it would be interesting to investigate the magnetization process of 
the S = 1/2 Heisenberg model on the Kagome lattice, where attention has so far been 
concentrated on the low-lying excitation spectrum (see [51] and references therein). 

In a final part, we have numerically computed the asymptotic behaviour of the magneti- 
zation curve close to the transition to saturation. The fundamental excitations associated 
to this transition are hard-core bosons and one therefore finds a universal behaviour: The 
characteristic DN-PT square-root in one dimension [38, 39], a linear behaviour with a log- 
arithmic correction in two dimensions [19, 35, 36] and a simple linear behaviour in three 
and more dimensions [42]. This issue was studied for the transition to saturation because 
the diagonalization simplifies considerably in this limit. In particular, one can explicitly 
map the problem to a low-density gas of hard-core bosons. In one dimension, where the 
nature of the excitations is not really important, almost all second-order transitions at 
plateau-boundaries are in the DN-PT universality class. So, presumably the universality 
class observed at the transition to saturation is also more general in higher dimensions. 
However, there is also room for different behaviour since now the nature of the fundamental 
excitations and the interactions are more important. 
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